import numpy as np  #导入numpy包


def EWMA(values, precision=1.e-3):             #def用于定义函数， 函数名（参数1，参数2，…）
    # 由给定精确度确定划分格点数。
    M = int(1 / precision)               #在[0,1]区间内均匀地取M个点，作为λ的尝试取值
    values = np.array(values)             #列表
    # 计算 {u_i}
    U = (values[1:] - values[:-1]) / values[:-1]    #计算日收益率
    U_squared = U * U

    opt_lbd = None
    min_loss = float("inf")      #正无穷
    # 穷举找最优lambda
    for i in range(1, M):
        lbd = float(i) / M            #λ=i/M
        sigma_squared = U_squared[0]   #设初始日的方差=u1²
        loss = 0
        for j in range(1, len(U_squared)):
            loss += np.log(sigma_squared) + U_squared[j] / sigma_squared    #目标函数
            sigma_squared = lbd * sigma_squared + (1 - lbd) * U_squared[j]    #EWMA模型的公式
        if loss < min_loss:
            min_loss = loss
            opt_lbd = lbd

    # 用最优lambda再计算出日方差率估计值。
    Vars = [0, U_squared[0]]
    for i in range(1, len(U_squared)):
        Vars.append(Vars[-1] * opt_lbd + (1 - opt_lbd) * U_squared[i])

    return (Vars, opt_lbd)

data = np.genfromtxt("Index50.txt")
Vars, lbd = EWMA(data, precision=1.e-3)
print("EWMA最优lambda：", lbd)
print(Vars[-1])
